Fermionic Molecular Dynamics 



H. Feldmeier* and J. Schnack^ 

Gesellschaft fiir Schwerionenforschung mbH, 
Postfach 110 552, D-64220 Darmstadt 



Abstract 

A quantum molecular model for fermions is investigated which works with an- 
tisymmetrized many-body states composed of localized single-particle wave 
packets. The application to the description of atomic nuclei and collisions 
between them shows that the model is capable to address a rich variety of 
observed phenomena. Among them are shell effects, cluster structure and 
intrinsic deformation in ground states of nuclei as well as fusion, incomplete 
fusion, dissipative binary collisions and multifragmentation in reactions de- 
pending on impact parameter and beam energy. Thermodynamic properties 
studied with long time simulations proof that the model obeys Fermi-Dirac 
statistics and time averaging is equivalent to ensemble averaging. A first order 
liquid-gas phase transition is observed at a boiling temperature of T ~ 5 MeV 
for finite nuclei of mass 16 ... 40. 
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1 Introduction and summary 



Heavy-ion reactions show typical dissipative phenomena at low beam energies (a 
few MeV per nucleon above the Coulomb barrier) |], §. With increasing impact 
parameter the complete fusion reactions go over into deeply inelastic collisions where 
the scattered nuclei have lost large fractions of the initial energy and converted that 
into intrinsic excitation energy. The exchange of nucleons causes the mass and charge 
numbers of the outgoing nuclei to fluctuate around a mean value. At the same time 
energy loss and scattering angle show fluctuations like a Brownian movement. For 
the description of these phenomena the so called particle exchange picture [[| |], |5[ 
which is based on the assumption of a Fermi gas with long mean free paths turned 
out to be rather successful. 

Around the same time the time-dependent Hartree-Fock model (TDHF) [[], [7], 
[|] was conceptually and numerically developed to a stage where heavy-ion collisions 
could be calculated. Although TDHF is a microscopic quantal model, which is 
supposed to describe systems with long mean free path and slow collective motion of 
the mean field, slow compared to the Fermi velocity such that the particles and their 
mean field can always be in equilibrium, it turned out that only the dissipation of 
the collective energy could be described but not the fluctuations || which inevitably 
go along with any dissipation. TDHF results are very close to a classical trajectory 
picture with friction and very small fluctuations. 

The surprising failure of the more quantal and more microscopic TDHF model 
compared to the more phenomenological particle exchange picture, both being based 
mainly on independent particle motion, originates in the inability of the TDHF state 
to react to small fluctuations. Consider a symmetric system: if there is a fluctuation 
in one direction, let us say a nucleon passes from the left to the right hand nucleus, 
there is of course in a many-body quantum state, always with the same probability 
the fluctuation in the opposite direction. Thus the mean field which averages over 
all configurations, will not adjust to the new situation, namely that there is now one 
more particle in the right hand nucleus, because with exactly the same probability 
a particle went out to the other nucleus and thus the mean particle number stays 
the same. The same holds true for recoil effects which would give a fluctuation 
to the collective relative momentum. The time evolution with the Hartree-Fock 
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Hamiltonian conserves global symmetries of the initial state. The mean field in 
TDHF is the same for all macro-channels (mass numbers, charge numbers, relative 



momenta etc.) which causes a "spurious cross channel coupling" [10| suppressing 
almost all fluctuations. 

In the particle exchange picture a random transfer of a nucleon will cause a jump 
to a new mean field situation where the acceptor nucleus is now enlarged in volume 
such that the new radius is that of a nucleus with one more particle. This way 
the decision is made, which of the many possibilities in the many-body quantum 
state is realized and this new situation is then evolved further in time. Of course 
with the same probability the particle could have jumped to the other side, then 
this would be the new situation with a larger mean field on the other side which 
would be evolved further in time. This picture is much more realistic for dissipative 
phenomena because the exact solution of the Schrodinger equation develops random 
phases between the macro-channels such that they cannot interfere on a macroscopic 
scale. 

When the beam velocity becomes comparable to the Fermi velocity the mean 
field picture of binary collisions should break down. In addition the mean free 
path of the nucleons becomes shorter than the diameter of the system due to the 
larger intrinsic excitation energy. In BUU-type models [II, |12|, [T^, [TJ|, [Tj| [RJ 



a. 



random Boltzmann collision term has been introduced to account for the shorter 
mean free path, but the long range part of the two-body interaction is still treated 
in a mean field fashion. The next step has been to replace the mean-field part 
by classical molecular dynamics with soft phenomenological two-body potentials 
fT/j , |j~8| . These Quantum Molecular Dynamics (QMD) approaches have two main 
conceptional difficulties. First, the mean field (long range part of the interaction) and 
the collision term (short range repulsion) should be treated self-consistently |H| |2(| . 
Second, the Pauli principle is reduced to the numerically enforced prohibition of 
over-occupation of phase-space cells. Other approaches try to replace the Pauli 
principle by adding a two-body "Pauli potential" in order to avoid too large one- 



body phase-space densities |2T], [22|, |3|, |J, £5], p6fl . The difficulty is that the Pauli 
principle requires the many-body state to be antisymmetric with respect to particle 
exchange, and this is outside the classical notions. 

In the following sections we discuss a model, named Fermionic Molecular Dy- 
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namics (FMD) |27], |29|, which combines the microscopic quantal features includ- 
ing the Pauli principle with the properties of classical molecular dynamics. Wave 
packets replace the classical points in phase space and thus introduce quantum prop- 
erties into classical molecular dynamics. Especially the antisymmetrization of the 
many-body state has many consequences genuine in quantum and absent in classical 
mechanics. 

A classical system of particles (molecules) which interact via a two-body inter- 
action of van der Waals type — repulsive at short and attractive at large distance 
- often behaves in a deterministic chaotic way. Small deviations in the initial con- 
ditions lead to exponentially diverging trajectories in the many-body phase space. 
Situations with global symmetries, for instance the left-right symmetry in the TDHF 
solution discussed above, are of measure zero in the classical case. Each particle fol- 
lows a single trajectory which cannot split into two or more with certain amplitudes 
like in the quantum case. This property is also common to FMD since there is only 
one wave packet per nucleon for which the equations of motion decide in a deter- 
ministic way where to move, for example whether it will join the right or left hand 
nucleus. This "quantization" of the particle number density is the main difference 
to TDHF if a single Slater determinant is used. In a forthcoming publication we are 
introducing antisymmetric many-body states which are not single Slater determi- 
nants anymore, but take care of the short range correlations caused by the repulsive 
core in the interaction. 

After discussing the concept of the model and the resulting properties we calcu- 
late ground states of nuclei. In section ^ an important insight is that the antisym- 
metrization not only takes care of the Pauli principle, but also delocalizes the wave 
packets and introduces shell effects. The ground states turn out to be rather similar 
to deformed Hartree-Fock states. 

The dynamical properties of the FMD equations of motion are investigated in 
section Calculations of heavy-ion collisions at 6 AMeV show that the typical 
fluctuations of dissipative collisions, which are absent in TDHF, are seen in FMD. 
Increasing the beam energy to 32 AMeV which corresponds to a relative velocity 
of about the Fermi velocity leads to a break down of the mean field. The very 
same model which describes fusion and binary dissipative collisions now predicts 
mult if r agment at ion . 
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A microscopic transport model developed for non-equilibrium situations should 
of course also possess the correct equilibrium properties. In the last section we 
proof that the thermostatic properties calculated with the FMD trial state are those 
of Fermi-Dirac statistics. Time averaging over dynamical calculations defining the 
thermodynamic properties not only shows that FMD equilibrates towards the Fermi 
distribution, but even allows to investigate the liquid-gas phase transition of small 
nuclei. 



2 The concept of Fermionic Molecular Dynamics 

The concept of Fermionic Molecular Dynamics is based on the molecular dynamics 
picture in which subgroups of the investigated system are described by their centre 
of mass coordinates. These subgroups are composite bound objects like molecules 
or in our case nucleons. The molecules are interacting by two-body potentials which 
approximate the complicated interactions between the constituents of the molecules. 
The two-body potentials, which include to a certain degree also the polarization 
induced in one molecule by the presence of the other one, depend in general on the 
relative distance, the relative velocity and the relative orientation of the molecules. 

Consider for example the interaction V ww between water molecules which have 
a dipole moment. First, it is a function of the distance but also of the relative 
orientation (see fig. p. Second, depending on the relative velocity the electron 
clouds will be able to adjust more or less completely to the adiabatic situation. This 
will result in a velocity dependence. 

In our case the complex molecule is the nucleon which consists of three quarks 
surrounded by a cloud of virtual mesons. Like the water molecule the nucleon is 
not spherically symmetric but has an orientation, the spin. In addition there are 
two kinds of molecules, protons and neutrons, which introduces a further degree of 
freedom, the isospin. Altogether the picture of a nucleon-nucleon potential in itself 
is only an approximation and there is a certain degree of freedom in its choice. 

In classical molecular dynamics the molecules are treated as classical distinguish- 
able particles which means that one has to solve Newton's equations of motion for 
the classical centre of mass coordinates. The many-body state is given by the tra- 
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Figure 1: Interactions between water molecules and nucleons (protons or neu- 
trons) are rather complex. 

jectory in many-body phase space. 

Here Fermionic Molecular Dynamics differs substantially because its many-body 
state is described by an antisymmetrized many-body wave function for the centre 
of mass coordinates and the spin degrees of freedom of all nucleons. As already 
explained in the introduction, in nuclei the nucleons cannot be localized in phase 
space well enough for a non-quantal treatment. The same holds true for the spin 
variables. Unlike the water molecule, where the dipole moment has only small 
quantum fluctuations around its mean value, the Hilbert space for the spin is only 
two-dimensional and hence one cannot construct a localized wave-packet for the 
spin orientation. A spin | degree of freedom is always quantal. 

All this together does not allow to treat nucleons in nuclei as classical particles. 
Therefore FMD is using localized single-particle wave packets for the centre of mass 
coordinates instead of points in phase space. The indistinguishability of the particles, 
which has such prominent consequences like Fermi motion, Fermi pressure, or Pauli 
blocking, is introduced by using a many-body state which is antisymmetric with 
respect to particle exchange. Therefore, all effects which arise from the fact that the 
particles are indistinguishable fermions are taken into account and by construction 
the Pauli principle cannot be violated. The use of wave packets guarantees that 
the uncertainty relation cannot be violated either. It actually turns out that due 
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to antisymmetrization the localization of the particles is suspended when the wave 
packets are overlapping. Thus, shell model states as well as plane waves are included 
in \Q) even if it is built from Gaussian shaped single-particle states (see section 



This means that the trial state contains the essential quantum properties of a 
Fermi gas. In addition if the system is dilute and the distance between the particles 
is much larger than the widths of their wave packets one obtains the classical limit 
in the sense of Ehrenfest. Hence the FMD trial state has enough freedom to cover 
large areas of quantum physics and connects seamlessly to classical physics. 

So the first part of the FMD concept is to include all necessary physics in the 
many-body trial state | Q ) but keeping this many-body state still numerically 
treatable. 

The second part of the concept is to deduce the equations of motion from the 
time-dependent quantum variational principle 

5 j\t{Q{t)\ ij t -H \Q(t)) = 0. (1) 

This variational principle fl30 |, which will be explained in more detail below, ensures 



automatically the conservation laws, provided the trial state has enough freedom as 



shown in section 2.3 



Altogether, the concept is well defined and the success of the model is granted if 
one has the proper Hamiltonian H and a rich enough trial state. The art consists in 
finding a trial state which contains the essential degrees of freedom without becoming 
numerically intractable. As one does not know a priori what the important physical 
phenomena will be, one needs experimental guidance and physical intuition to set 
up the trial state. Even for the same Hamiltonian H there can be different optimal 
sets of trial states, depending on the excitation energy or the part of the phase space 
the system is occupying. 

If the model fails, the concept is to reexamine the trial state and to include 
more or other degrees of freedom. But also the Hamiltonian is not sacred as it is in 
many cases only an approximation to a much more complex underlying microscopic 
picture. 

There is also an interplay between the trial state and the effective Hamiltonian. 
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In nuclear physics the nucleon-nucleon interaction is strongly repulsive at short 
distances. This important feature contributes to nuclear saturation and hence is 
part of the essential physics which should also be present in the trial state. The 
antisymmetrized product state (single Slater determinant) which we use in this paper 
cannot account for the depletion of the many-body wave function at small relative 
distances between the particles. Therefore, we are using an effective interaction 
with a moderate repulsion at short distances in the spirit of a Brueckner G-Matrix 
treatment. 



2.1 Parameterized many— body trial state 

Instead of classical points in phase space Fermionic Molecular Dynamics deals with 
wave packets which are localized in phase space. Each nucleon (molecule) is repre- 
sented by a Gaussian wave packet parameterized by the set 
q(t) = {7>(t),p(t),a(t),x(t),<j>(t),Z}: 

(x \q(t)) = exp|- 2fl( y ; +ip{t)-xj® \x(t),(f>(t))® |£) , (2) 

which in quantum mechanics is the closest analogue to a classical particle described 
by a single point (f(t),p(t)) in phase-space. 

In FMD the match to classical positions and momenta are the parameters f(t) 
and pit) which determine the mean values of the position and momentum operator 
of the single particle state: 

(g(t)\ jU(i)) (sWljkW) 

r(t) = mum • m = (W> • (3) 

Due to the quantum mechanical uncertainty relation the wave packet can be either 
narrow in coordinate space and wide in momentum space or vice versa. This non- 
classical degree of freedom is taken care of by the complex width parameter a(t) = 
a R,(t) + ia>i(t). It determines via its real part a#(t) the variance of the momentum 
distribution a 2 K {t) by the relation 

(l(t)\(k-p(t)) 2 \<l(t)) „ 2 



2a R (t) (q(t)\q(t) 



= 3<(t) • (4) 
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Since the wave packet is spherical the widths are equal in all three Cartesian di- 
rections. The imaginary part ai(t) appears in the expression for the spatial width 
(Tx(t) as 

3 al(t) + a]{t) _ ( q{t) | (g - r{t)? I g(*) ) _ 

2 a*(t) < 9(0 I 9(0 > X 

and determines in how far the wave packet is of minimal uncertainty. The product 
of the variances 

«4 4 = \ (l + }) (6) 

shows that for a/ = one has the minimum-uncertainty packet where <Jx<Jk = \ , 
while for ai ^ the uncertainty can become arbitrarily large. This means that the 
particle occupies more than (h/2n) 3 of phase space volume but at a lower phase 
space density, such that other fermions can find place at the same area in phase 
space. 

Besides the parameters for the spatial part of the wave packet there are two 
parameters x(t) and (pit) for the spin degree of freedom. If one parameterizes the 
trial spin-state by 

<™.ix(*),w)>=( cos x f m) : ms= K . ( 7 ) 

{ sin^eW : m s = -\ 
the relation between the parameters and the corresponding spin operators is 

(q(t) | a \q(t)) 

where a(t) = (sin%(i) cos0(t), sinx(t) sin0(t), cosx(t)) is a vector in the 3-dimensional 
real space and the quantization axis is the 3-axis. Thus a(t) can be regarded as 
the "classical" spin direction just like f(t) or p(t), although there is no classical spin 
degree of freedom which can vary only its direction but not its magnitude as it is 
the case for a(t). 

In principle the same parameterization can be chosen for the isospin part | £ ). A 
time-dependent isospin would mean that for example due to the exchange of charged 
pions neutrons can dynamically transform into protons and vice versa. Isospin sym- 
metry is, however, only an approximate symmetry, the Coulomb interaction and the 
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difference in proton and neutron mass break it. A linear superposition of protons 
and neutrons in one wave packet would impose undesired artificial symmetries. For 
example only the proton component would feel the Coulomb repulsion and its cen- 
tre of the wave packet would be accelerated away from the centre for the neutron 
component. Therefore, either one gives each component its own wave packet or one 
does not allow proton neutron mixing in the trial state. In this paper we shall not 
consider rotations in isospin space but assume | £ ) to be independent of time and 
either | £ ) = | proton ) or | £ ) = | neutron ) . 

To describe a system with A fermions we construct a Slater determinant | Q(t) ) 
with these parameterized single-particle states | quit) ) 

\ qi (t),q 2 (t),--- ,q A (t))= |Q(t)) = _^_L_ |Q(I)), (9) 

(Q(t)\Q(t))* 

where the antisymmetrized but not normalized state | Q{t) ) is given by 

all tt 

The sum runs over all permutations it and sgn(7r) is the sign of the permutation. It 
should be noted that quit) = {rk{t),Pk{t),<ik{t),Xk{t),<pk{t),£k} denotes the set of 
parameters specifying the single-particle state with number k. The parameter set 
for the many-body state thus reads 

Q(t) = {nit),^),^^)^!^),^)^!,^),--- ;r A (t),--- ,£ A } 

= {q v {t)\v = l,--- ,NA} , (11) 

where N is the number of real parameters per particle, in our case N=10. Whenever 
q carries a Greek index it is an individual parameter, whereas a Latin index implies 
that qk is the whole set for the state \qt)- 

Due to antisymmetrization FMD is constrained to the antisymmetric subspace 
of the Hilbert space and hence the Pauli principle is a priori incorporated. Further- 
more, the projection (|1(]) from a product state onto the antisymmetric subspace 
destroys for overlapping Gaussians the localization of the particles and introduces 
shell model states. This will be discussed and explicitly shown in section |3.2| . If 
the single-particle states | qk ) are not overlapping, the antisymmetrization has no 
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effect anymore and the particles are localized in the individual wave packets. In this 
limit we return to classical Newtonian mechanics for f&(i) and Pk(t), which however 
can still be coupled to the non-classical variables <7fc(t) for the spin directions and 
the widths dk(t). Of course also in this limit the particles are indistinguishable and 
it is not possible to decide which particle occupies which Gaussian packet. 

Even though eqs. (|3|) - fl5|) and (§) provide unique relations between the param- 
eters and the expectation values of the corresponding operators, this is no longer 
true for the antisymmetrized many-body state. To illustrate this important aspect 
let us regard a two-body state of two packets which for simplicity have the same 

real width parameter a\ = a 2 = a# and the same spin, but different f*i, r 2 and pi,p 2 

I 

I 91,92)0 = I 9i) ® I 92) - I 92) ® |9i)} • (12) 



/2N 

with the normalization 

^=(9i|9i>(92|92>-|(9i|92>| 2 (13) 
= (91I91K92I92) (l-exp{-f? 2 }) , ^ = (fx -f 2 ) 2 /a R - (pi -p 2 ) 2 a R , 

where £i 2 measures the distance in phase space. 

The operator for the position of particle 1, ^(1) = % <8> 1 , (which by the way is 
not an observable because it is not symmetric under particle permutation) has the 
following expectation value 

a( 91, 92 I £(1) |9l,92)o (14) 

= 2^{<9il i I 9i ) + (92 I 1 I 92) -2Re((g x | x | 92 ) ( 92 I 9i ))} 
1 -* N 

= j ( r l + r 2) • 

The result is not ri as presumed from Ehrenfest type arguments, however, 
(£(1) + £(2)) = n + r 2 . 

Another example is the relative distance between two identical fermions, x(l) — 
x (2), (again not an observable, but used for example in its classical meaning in the 
collision term of QMD ||18||) with the expectation value 



a(9i,9 2 | £(l)-£(2) \ qi ,q 2 )a = 0. (15) 

If the spin of the two fermions is the same, the two-body state has to be antisymmet- 
ric in coordinate space and hence this expectation value is always zero, independent 
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on the actual value of T\ and f 2 . The same holds true if the spins are anti-parallel. 
If the spins are not parallel the "distance" depends on the spin orientations, which 
only reflects that f\ — f 2 is not the mean distance between the two fermions. 

An observable is for example the rms-distance 1^(1) — x (2)| which, however, 
also depends on the spin directions and the relative momentum. The reason is 
that for the component of the trial state in which the two fermions have total spin 
S — the radial part of the wave function is symmetric, while for the three S = 1 
components it is antisymmetric. Hence the rms-distance is different for the singlet 
and the triplet component. In addition the exchange term contributes less if the 
relative momentum becomes large. The Pauli principle, which demands that two 
fermions cannot be in the same phase space cell, is of course the origin of these 
classically not existing correlations in the rms-distance. 

Another example for an observable which is even used to motivate a so called 



"Pauli potential" is the kinetic energy [p3| , [24]]. For equal spins and equal widths 
oi = &i = clr we get 

Z)~a(quq2\k\l) + k 2 (2)\q h q 2 ) a (16) 

I l r -2 -2 -2 ^ 

= 2^ivV gi| ~ \to) + (Q2\k -2Re(( 9l | k | q 2 )( q 2 | 9l )) j 

1 | {Pi ~ Vif + {n - f 2 f/a R exp {-^ 2 2 } 



2m l-exp{-^ 2 2 } 



Am 2maR 

One sees that the kinetic energy is not simply (pf +pf)/(2m) but depends also on 
the "relative distance" fx — f 2 and on the distance in phase space £12. But one has 
to be aware that for wave packets which are not of minimal uncertainty (i.e. the 
imaginary part aj of the width parameter a is not zero) the whole expression is 
completely different and |( q 2 \ q\ )| 2 can be zero even when the packets fully overlap 
spatially. A non-vanishing a/ is nothing exceptional, it describes for example the 
well-known spreading of the packets in free space. Therefore, the constraint aj = 
is not appropriate for the dynamical case and leads actually to unphysical scatterings 
and strong hindrance of evaporation. 



In section |3.2| we shall show for an A-body system that even if all parameters for 
the momenta vanish, i.e. p\ — p 2 — • ■ ■ — pa — 0, the kinetic energy is neither zero 
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nor the sum of the zero-point energies 3/(2ma#), but the momentum distribution 
is that of a Fermi gas with a sharp edge at the Fermi momentum. 

To summarize this discussion, one should always carefully distinguish between 
parameters and physical observables. The parameters are in general not the expec- 
tation values of the corresponding operators. 



2.2 Time— dependent variational principle 

After having set up the many-body trial state \Q(t)) one has to construct equations 
of motion for the set of parameters Q(t) = { q u (t) \ v = 1,2,3,- •• } which are 
the generalized coordinates of the system. As shown in the previous section the 
parameters may loose their original physical meaning due to antisymmetrization. 
Therefore one must not simply use the classical equations of motion for the position 
and momenta. In addition, even without antisymmetrization, classical physics does 
not tell how the complex width parameter should be evolved in time. Also the 
equations for the angles (xk, 4>k) of the spin direction a, which turn out to be of the 
Bargmann-Michel-Telegdi type fl3l |, are not self-evident. 



Therefore, the equations of motion for the Fermionic Molecular Dynamics model 
are derived from the following time-dependent variational principle 

5 j\(Q(t)\ ij t -H \Q(t)) = (17) 

in which the trial state | Q(t) ) is to be varied. This has the advantage that the 
geometry of the manifold of trial states is automatically taken care of, irrespective 
of how the parameters of the states are defined. Furthermore, conservation laws 
follow in a transparent way from invariance properties of the trial state and the 
Hamiltonian, see section ET3". 



The variation has to be performed with respect to each parameter q u {t) with 
the end points kept fixed, i.e. 8q v {t\) = Sq u (t 2 ) = 0. The operator H is the total 
Hamiltonian of the system. (Throughout the paper operators in Hilbert space will 
be underlined with a twiddle to distinguish them from parameters or expectation 
values.) 
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The Euler-Lagrange equations 
d dC dC 



, v = l,2,---,N (18) 



dt dq u dq v 

which result from the variation flT7| ) are written in terms of the Lagrange function 

c(Q(t),Q(t)) := (Q(t)\ ij t -H \Q{t)) 

= C (Q(t),Q(t)) - H(Q(t)) , (19) 

with 

Co (Q(t), Q(t)) := ( Q(t) | i- | Q(t) > = ^ ( Q(t) 1 1— | Q(t) ) q u , (20) 

in which Q(t) = {q u {t) = dq v /dt \ v = 1,2, ••■} is the set of generalized veloci- 
ties and 7i(Q(t)) is the Hamilton function defined as the expectation value of the 
Hamiltonian H\ 

H(Q(t)) = (Q(t)\H\Q(t)) . (21) 

Different from classical mechanics the Lagrange function (^) is linear in the veloc- 
ities q u but at the same time the set Q(t) = {q u {t)} contains both, coordinates and 
momenta. 

Using the general structure ( |T9"D of the Lagrange function the Euler-Lagrange 
equations ( |I8"D in their most general form can be written as 

EA.(«)i = -^ (22) 



or, if A^v is not singular |!2| the equations of motion are 



?» = -£^(Q)^P. (23) 

V 

where 

= -A,(Q) = - , (24) 

is a skew symmetric matrix, which depends in general on all variables Q(t) = {q u (t)}. 
For details see ref. p8 . 
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With help of the matrix A^ v one can define generalized Poisson brackets f3(| as 

such that the time derivative of an expectation value 

B(t) = (Q(t)\B\Q(t)) (26) 
of a time-independent operator B calculated with the trial state | Q(t) ) is given by 

= Eg^'f = <«■">■ (*) 



Equation fl27|) has the symplectic manifold structure of Hamiltonian dynamics |33| 



but in the general case the parameters q u cannot be grouped into pairs of canon- 



ical variables. However, according to Darboux's theorem [33] canonical variables 
exist locally. They are non-linear functions of the parameters q u and have to be 
constructed such that A~l assumes the canonical form 

where 1 is the unit matrix. Their choice, however, is not unique. 



2.3 Conservation laws 

After having solved the equations of motion (^) for the parameters Q(t) = {q u {t) \ v - 
1, ■ • • , N} the trial state | Q(t) ) is known at all times. Thus one can calculate the 
expectation value Q(t) = (Q(t)\ G \Q(t)) of an arbitrary time-independent opera- 
tor G. With the definition ( p5]) for the Poisson brackets the time derivative of this 
expectation value can be written as 

J t S = {H,G}. (29) 



The expectation value is conserved in time if ]33 
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Since A^l is skew symmetric the energy 7i itself is always conserved by the equations 



of motion, provided they are derived from the variational principle fll7|) . This is 
completely independent on the choice of the trial state . 

For other constants of motion we consider a unitary transformation with G as 
the hermitean generator 



U = expyieG\ , e real. (31) 
If U maps the set of trial states onto itself 

U\Q)= \Q') e {\Q)} , (32) 



then, as a result of the equations of motion, it can be shown [28[ that 



{H,g} = (Q(t)\i[H,G\ \Q(t)) . (33) 

This means that for this class of generators the generalized Poisson bracket is just 
the expectation value of the commutator with %H. 

Relation (^) is very useful for two reasons. First, if G commutes with the 
Hamiltonian H and exp(ieG) | Q ) = \Q') then ( Q(t) \ G | Q(t) ) is automatically a 
constant of motion. 

Second, this relation is an important guidance for the choice of the trial state 
\Q). If one wants the model to obey certain conservation laws then the set of 
trial states should be invariant under the unitary transformations generated by the 
constants of motion. For example, total momentum conservation implies that a 
translated trial state is again a valid trial state. This is fulfilled for the trial states 
specified in section |27I|. The Gaussians defined in eq. (0) can be translated or Galilei 



boosted, the latter taking care of the conservation of the centre of gravity. 

Conservation of total spin J = L + S is guaranteed if rotation of the trial state 
in coordinate and spin space results again in a trial state. This implies that the 
Gaussian (^) has to have either a spherical shape or the width parameter has to be 
replaced by a complete tensor with 12 real parameters. It also means that all spin 
directions in \4>,x) have to be allowed, otherwise the rotation would in general lead 
out of the manifold of trial states. 

If G does not commute with H : relation ( |3~3D sheds some light on the quality of 
the variational principle (|17D . It says that under the premise that exp(zeG) does not 
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map out of the set of trial states the time derivative of the expectation value of G 
calculated with a trial state is exact. 



* 6(t) = 5 6 — (i) • (34) 

where the exact solution with the initial state | Q(t) ) is 

<3e Xac t{t + t) = ( Q(t) | e lT ~ G e~ lT ~ \ Q(t) > . (35) 

Due to the fact that | Q(t + r) ) is only the approximate time evolution of | Q(t) ) 
for some time r later Qit + r) will begin to deviate from Q e xact{t + t) for larger r. 

The kinetic energy is such a generator. Since our trial state (§) is the exact 
solution of the Schrodinger equation without interactions it fulfills 

eM-irT)\Q(t))= \Q(t + r)). (36) 

With the two-body interaction included, | Q(t) ) is not an exact solution anymore, 
but the expectation value of the total kinetic energy, which then is not a conserved 
quantity any longer, is well approximated in the sense of equation (|34|). 



2.4 Two— body Hamiltonian 

2.4.1 Effective nucleon— nucleon interaction 

Up to now the interaction V contained in the Hamiltonian was not specified. Since 
in this paper we shall investigate small and medium-heavy nuclei, we choose an 
effective two-body potential suited for mass numbers up to about 50. This inter- 
action is repulsive at small and attractive at larger distances. The repulsive core 
is, however, rather weak and one should regard the potential as a phenomenological 
ansatz for a G-matrix rather than the free nucleon-nucleon interaction which has a 
very strong repulsion for distances smaller than 0.5 fm. In this paper we are using 
a potential of the form 

V(i,j) = V a (w a + m a P r + b a P° + h a PT) (37) 
+ Z b K + m b P r + b b P° + h b P T ) 
+ V c ( Wc + mc Pj + bc P° + h c Pj) 
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where P r ^P (T ^Pj denote the various exchange operators for coordinate, spin and 
isospin. The radial dependences are of Gaussian type given by 

I Z h \%k,xi) = V aAc exp \ - fe 2 X ^ \ (38) 

{ r *Ac J 

x5 3 (xi - x k )5 3 (xj - f/) 

which allows to calculate all matrix elements and their derivatives analytically. The 
parameters 

K = -7.10MeV; r a = 1.16fm (39) 
V b = -31.90MeV ; r b = 2.22fm 
V c = +81.65MeV ; r c = 0.735fm 

w a = 8.700, m a = 5.610, b a = 7.860, h a = -21.170 

w b = 0.133, m h = 0.514, 6 fe = 0.085, h b = 0.268 

w c = 1.000, m c = 0.000, b c = 0.000, h c = 0.000 

were determined in order to reproduce the binding energies of a wide range of 
medium-heavy isotopes. 



2.4.2 Coulomb interaction 

The Coulomb interaction is included in the Hamilton operator, it is given by 

{x i ,x j \V{i,j)\x k ,x l ) = I^M^ P p PV(f, -x k )5\x, -£,) , (40) 

where P P denotes the projection operator on the protons 

P p \q m ) = \{l + U)\q m )- (41) 

The isospin variable £ m takes values £ m = 1 for protons and £ m = — 1 for neutrons, 
respectively. 

In order to speed up computing time the expectation values of the two-body 
potentials are approximated as given in the appendix and spin degrees have been 
kept fixed in the following calculations. 
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3 Ground states in FMD 



The ground state of a nucleus is the many-body state | Qgs ) in which the energy 
7~L — ( Qgs I H I Qgs ) is an absolute minimum with respect to variation of all 
parameters q v , therefore 

^^ = 0. (42) 

This implies that the FMD ground state is completely time-independent (up to 
an overall phase) and the time derivatives of all parameters vanish because, by 
definition, all generalized forces dTi/dq^ are zero and hence 

^ = -E<jg = °- ( 43 ) 

Requirement ( f4"2|) not only determines the positions fk and momenta f>k but also the 
complex widths and the spin directions (xk,<fik)- 

Since we do not correct for the centre of mass motion in the dynamical calculation 
we use the expectation value of the Hamilton operator as the ground state energy 
and do not subtract the centre of mass energy. The centre of mass energy is of the 
order of lOMeV for all isotopes, so that its contribution to the energy per particle 
vanishes for larger mass numbers. An ansatz for the many-body state where centre 
of mass and relative motion separate is introduced in ref. |35| . 



Besides the ground state energy we also calculate the root-mean-square radius 
of the charge distribution. 

E GS = (Qgs\H\Q gs ) (44) 

Rims = ^J2(Qos\&)-X c J^^\Q GS ) + Rl roton , (45) 

where the centre of mass position operator is 

A 



X 

i=l 



The operator P P (i) projects on protons and R pro ton = 0.876 fm takes the finite 
charge radius of the proton into account. 



Fermionic Molecular Dynamics 



21 



3.1 Ground states of nuclei 



In table [I] we summarize ground state binding energies and charge radii and compare 
them with experimental results. The interaction ( j37l) describes a wide range of 
isotopes with a satisfying accuracy. 





experiment 


FMD 


experiment 


FMD 


Isotope 


Egs 


Egs/A 


( H ) 




Rrms 


Rrms 




( MeV ) 


( MeV ) 


( MeV ) 


( MeV ) 


( fm ) 


( fm ) 


4 He 


-28.296 


-7.07 


-28.32 


-7.08 


1.63 


1.63 


12 C 


-92.163 


-7.68 


-91.36 


-7.61 


2.42 


2.79 


16Q 


-127.62 


-7.98 


-127.7 


-7.98 


2.73 


2.84 


19 F 


-147.80 


-7.78 


-143.3 


-7.54 




2.87 


24 Mg 


-198.256 


-8.26 


-196.7 


-8.19 


2.95 


3.17 


27 Al 


-224.952 


-8.33 


-216.9 


-8.03 


2.95 


3.13 


28 Si 


-236.537 


-8.44 


-238.3 


-8.51 


3.04 


3.24 


32g 


-271.783 


-8.49 


-273.4 


-8.54 


3.20 


3.35 


40 Ca 


-342.056 


-8.55 


-349.0 


-8.72 


3.50 


3.44 



Table 1: Ground state energies and charge radii in FMD. The experimental bind- 

and the charge radii from ref. |38|, p9|1 . 



ing energies are taken from ref. |3q , 
All experimental errors are at most in the last digit. 



In fig. and fig. [5] we display the density of different ground states in coordinate 
and momentum space, respectively. The densities are defined as 

p x (x) = (x\ P {1) \x) and p k (k) = (k\ P (1) \k) , (47) 

where pj 1 ^ is the one-body density operator. The crosses indicate the centres fi and 
Pi of the wave packets. One should however keep in mind that a Slater determinant 
is invariant under transformations among the occupied single-particle states. In 
the coordinate representation the density is integrated over the ^-direction, the 
momentum representation is shown as a cut at k z = 0. 

In fig. |2| one sees that the ground states are intrinsically deformed. The 12 C 
nucleus for instance arranges as three a-clusters, whereas this a-symmetry is broken 
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Figure 2: Contour plot of spatial densities integrated over z-direction for differ- 
ent ground states. Crosses indicate centroids of wave packets. Contour lines are 
at 0.9, 0.5, 0.1 and 0.01 of the maximum density. 
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Figure 3: Contour plot of momentum distributions cut at k z = for different 
ground states. Crosses indicate centroids of wave packets. Contour lines are at 
0.03, 0.06, 0.09, 0.12 and 0.15 fm 3 . 
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in 19 F and 27 Al. The true ground state, which is an eigenstate of the total spin, is 
a superposition of all orientations of the intrinsically deformed ground state |28] . 

The momentum distributions shown in fig. [3] reflect the Fermi motion of the 
ground states. Although all single-particle wave-packets are stationary (and even 
centred at p — 0) the system possesses Fermi motion. In FMD the Fermi mo- 
tion is a quantum mechanical zero-point motion and not a random motion of the 
packet centroids. In the ground state r k and Pk and all other parameters are time- 
independent otherwise it would be not the ground state of the system. Using a 



product state for the many-body system |L7|, [18], |20|, [40], [41], [42| (43[ it might be 
possible to use a momentum dependent Pauli potential in order to mock up the 
momentum distribution. Nevertheless such a system has the thermodynamic prop- 
erties of distinguishable particles, which are different from those of a Fermi system 
especially at low excitations. 



3.2 Shell structure in FMD 

It is not immediately obvious that FMD includes shell-model features like the nodal 
structure of single-particle orbits since the states are localized in coordinate and 
momentum space. But due to the invariance of a Slater determinant under linearly 
independent transformations among the occupied single-particle states, after anti- 
symmetrization, any set of single-particle states which is complete in the occupied 
phase space is as good as any other. This applies also to non-orthogonal states. 
To illustrate this we take four one-dimensional real Gaussians with the same real 
width parameter a and zero mean momentum and displace them by d = 0.75^ s /a 
(see l.h.s. of fig. ^). The one-body density can be written in terms of orthonormal 
states | ip m ) as 

A A 

P (1) =Y1 \lk)O kl ( qi \ =]T|^ m >(^| , (48) 

k,l=l m=l 

where the orthonormal eigenstates of P^ 1 ' are given by 

A 

\i>m) = J2 \9k)(0^) km (49) 

k=l 

and Oki is the inverse of the overlap matrix (qk\qi)- They are displayed on the right 
hand side of fig. f| and compared to harmonic oscillator eigenstates (dashed lines). 
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Figure 4: Antisymmetrization of the four Gaussians on the left hand side leads to 
harmonic oscillator states. Dashed lines are the exact eigenstates of the oscillator. 

One sees that the occupied single-particle states | ipm ) consist of an s-, p- ,d- and 
an f-state, all very close to harmonic oscillator states. The difference between both 
sets can be made arbitrarily small by letting d/y/a approach zero. 

A second example is illustrated in fig. [5], where we consider 100 equally spaced 



Gaussians in one dimension ||44|| . Again all mean momenta are zero and the width 
a is real. In the upper part of fig. [5] the width \fa is 0.2 of the mean distance d 
so that the wave packets are well separated. Therefore the spatial density p x and 
the momentum density p& are not changed by antisymmetrization. In the lower 
part the width has been increased to y/a = d. Without antisymmetrization (dash 
dotted line) the spatial density is uniform and the momentum distribution is that 
of a single packet. After antisymmetrization (full lines) one obtains the typical 
shell model oscillations in coordinate space and a Fermi distribution in momentum 
space. It is amazing to see how in eq. (|48|) the superposition of Gaussians by means 
of the inverse overlap matrices can create a fully occupied momentum state, see for 
example in fig. |5] the lower right momentum distribution at k — 0.8kp, where the 
individual Gaussians give practically zero probability to measure this momentum. 
We also calculated the eigenstates of the kinetic energy in the occupied space and 
got perfect sinusoidal waves. 

These two examples illustrate nicely that even localized single-particle states 
with zero mean momentum build up FMD many-body trial states which describe 
the harmonic oscillator shell model or even the Fermi motion of a gas of fermions in 
which plane waves are occupied up to the Fermi momentum. 
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Figure 5: Upper part: section of spatial density of hundred Gaussians (not over- 
lapping in coordinate space) and corresponding momentum distribution (same 
for all). Antisymmetrization does not change the distribution. Lower part: same 
as above but for overlapping Gaussians. Full line with antisymmetrization, dash- 
dotted line without. For details see text. 



If one wants to simulate this effect by a "Pauli potential", disregarding the 
momentum distribution in each wave packet, the resulting ground state momentum 
distribution is unsatisfactory p4]. 



4 Heavy ion reactions 

It is the aim of Fermionic Molecular Dynamics to describe the phenomena seen in 
heavy-ion collisions at laboratory energies of i?^,200 AMeV, as there are fusion, 
incomplete fusion, dissipative reactions, evaporation of nucleons and fragmentation. 
In order to gain a predictive power that goes beyond global phenomena dominated 
by conservation laws three basic conditions are indispensable. 
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The first condition concerns the ground state properties. Many observables like 
for instance fragment multiplicities depend strongly on ground state energies and 
radii. Thermodynamic properties like the specific heat are related to the exchange 
symmetry of the many-body state. Equilibration is strongly influenced by the mean 
free path of the constituents, which is usually much larger for Fermi systems than 
for classical systems. 

The second condition is that the trial state has to have the necessary degrees of 
freedom for the phenomena one wants to account for. Mean positions and momenta 
of the single-particle wave-packets are obvious degrees of freedom. For reasons of 
simplification the width is often chosen as real and time independent, but it could 
be shown that this degree of freedom is not only useful to reproduce the free motion 
exactly, but is needed to describe evaporation of nucleons and fragmentation of 
nuclei H, ||. 



The third important condition, which is related to the second, refers to sym- 
metries. The deterministic equations of motion of FMD preserve all symmetries 
in the initial state under which the Hamilton operator is invariant. It is therefore 
desirable to use a Hamilton operator that breaks as many symmetries as possible. 
This Hamilton operator should not only contain a central potential, as it is the case 
in the present calculations, but also spin-orbit and tensor interactions, which are 
expected to break spin symmetries of the ground states during the dynamics. In 
addition short range correlations have to be considered, which lead to more momen- 
tum transfer in the reaction. These short range correlations are responsible for hard 
collisions of nucleons, a role that is played by the fluctuating collision term in QMD 
or AMD. First steps into this direction are already taken p6 . 



The model is chosen to be a molecular dynamics model because it addresses 
the large fluctuations observed especially in multifragmentation reactions. In the 
model the molecular aspect is expressed through the localization of single-particle 
wave-packets; this can be regarded as a quantization of the particle number at any 
time. It means, fragments have always integer particle number, in contrast to time- 
dependent Hartree-Fock, see for instance f47| , f48fl . Another aspect is that the FMD 
single-particle wave-packets are not allowed to split, therefore tunneling through 
a barrier cannot be described. But this quantal process is of less importance in 
fragmentation reactions because it takes much longer than the involved time scales. 
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The equations of motion of FMD are deterministic, therefore the event ensemble 
does not arise from random fluctuations of the collision term during the time evolu- 
tion, but from the average over all orientations of the intrinsically deformed ground 
states. Although these orientations superimpose coherently in the true ground state 
their relative phases are randomized during the collision so that they add up inco- 
herently in the exit channels. The same holds true for the summation over different 



impact parameters [Eq, 49 



4.1 Deeply inelastic reactions 

As a first example for heavy-ion collisions described with FMD a dissipative reaction 
is presented. The reaction of 19 F and 27 Al at a laboratory energy of 5.9 AMeV was 
investigated at the SMP Tandem accelerator in Catania f50fl . At this low energy the 
reaction is dominated by two processes. For impact parameters b up to about 5 fm 
the system fuses completely, for larger b the reaction shows the typical dissipative 
phenomena of deeply inelastic collisions. The experiment focussed on the latter 
process. In order to understand the reaction a view on the nucleus-nucleus potential 
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Figure 6: 19 F- 27 A1 static nucleus-nucleus potential. The solid line shows the 
FMD result using the nucleon-nucleon interaction ([IT]) and frozen ground states, 
the crosses indicate the parameterization by Krappe, Nix and Sierk. 



might be helpful. In fig. [| the solid line displays the total energy of the 19 F- 27 A1 
system (less the ground state energies) as a function of the distance between the 
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two centres of mass as it arises from the two-body interaction fl37|) and (40|). The 
many-body state of the two nuclei was taken to be the respective ground states and 
not changed as a function of the distance. Therefore, antisymmetrization induces a 
strong repulsion below d — 5 fm. Above d = 5 fm this static interaction compares 
nicely to the parameterization of the nucleus-nucleus interaction by Krappe, Nix 
and Sierk |5l| , which was adjusted to bulk properties of nuclei and successfully 
applied in earlier studies of heavy-ion dynamics, e.g. in ref. @. 
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Figure 7: d 2 a/(d6 cm dTKE) for 19 F- 27 A1 collisions at 5.9 AMeV: FMD calcula- 
tions are displayed on the l.h.s., experimental results on the r.h.s.. Subsequent 
contours differ by a factor of 5. 



In the impact parameter range of 4.8 fm < b < 10 fm 509 events were generated 
with FMD. The results were filtered with the experimental angular acceptance 3° < 
#Lab < 54°. Below b = 4.8 fm the system was fusing. 

Figure shows the double differential cross section d 2 a/ (dd cm dTKE) as a function 
of the centre of mass scattering angle 8 cm and the total kinetic energy TKE of the 
two scattered nuclei. On the left hand side the FMD calculations are presented as 
a contour plot. Each event (# cm ,TKE) contributes as a Gaussian (A|# cm | = 2° and 
ATKE = 2 MeV) in order to smoothen the distribution. For impact parameters 
smaller than the one leading to a grazing collision (|# C m| ~ 15°, TKE pa 65 MeV) 
the approaching nuclei feel the nuclear attraction, 6 cm gets smaller and becomes 
negative. The total kinetic energy TKE is decreasing due to internal excitation of 
the nuclei. The major part of the cross section is observed at negative angles. For 
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even smaller impact parameters the nuclei stick together so long that they appear at 
positive angles again. One observes a flat distribution in 8 cm at the Viola energy of 
TKE « 20 MeV. On the right hand side the experimental result is given. Comparing 
both sides of fig. [7| one sees that FMD is capable to describe the dissipative character 
of the reaction dynamics. 
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Figure 8: Diffusion plot d 2 a/{dZ dTKE) for 19 F- 27 A1 collisions at 5.9 AMeV: 
FMD calculations are displayed on the l.h.s., experimental results on the r.h.s.. 
Subsequent contours differ by a factor of 5. 

The diffusion plot d 2 a/(dZ dTKE), fig. || samples fragments in the mentioned 
angular range of 3° < < 54°, which are mostly projectile like. For the experi- 
ment, r.h.s. of fig. [8], the data are limited by an energy threshold seen at low TKE 
and high Z, which is not imposed on the analysis of the FMD simulations. The 
FMD events (Z, TKE) are smoothened with A Z = 0.5 and ATKE = 2 MeV. 

Both contour plots, FMD (l.h.s.) and experiment (r.h.s.), show the typical broad- 
ening of the charge distribution with increasing dissipated energy and the lack of 
charge drift for TKE above the Viola energy of TKE « 20 MeV for the completely 
relaxed events. A drift to smaller Z values, which means that the 19 F nucleus is 
giving away nucleons to the 27 Al, is seen in the measured and the calculated cross 
sections for large energy losses where all the initial kinetic energy is dissipated and 
the reaction lasts very long. 

Figure |9] presents the second part of the diffusion plot for Z > 12 which has 
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Z 

Figure 9: Diffusion plot d 2 a/(dZ rfTKE) for 19 F- 27 A1 collisions at 5.9 AMeV 
calculated with FMD for target like fragments (Z > 12). Subsequent contours 
differ by a factor of 5. 

not been measured. Besides the deeply inelastic reactions one now also sees at TKE 
values close to zero the evaporation residues from fusion reactions with charges 
ranging from Z — 16 to the total charge of Z = 22 with a maximum around Z = 20. 
Here one should note that the FMD time evolution was followed only up to 1200 fm/c 
so that the loss of charges due to evaporation was not complete yet. Nevertheless 
the broad distribution of evaporation residues shows that FMD is not unrealistic in 
this respect. 

4.2 Multifragmentation 

In this section the same system of 19 F and 27 A1 is investigated, but at an energy of 
32v4MeV where multifragmentation is expected. 

How multifragmentation happens in heavy ion collisions is still a matter of de- 
bate. Explanations reach from nucleation over self organization, spinodal decom- 
position to cold break-up and survival of initial correlations. For an overview see 
ref . . A key question is the time scale of the reaction. Slow processes like nu- 
cleation or self organization are hindered if the expansion of the whole system is 
too fast. Another issue is the relaxation time for thermal and chemical equilibrium 



Fermionic Molecular Dynamics 



31 



which is important when statistical models are used to explain multifragmentation 



53, 54 



If one considers the decay of excited spectator matter which has not been com- 
pressed one expects a slow expansion so that there might be enough time for global 
equilibration. On the other hand the excitation energies are not so high, so that 
the mean free path, due to the Pauli principle, is still not small compared to the 
diameter of the nuclear system. Therefore it is not obvious that global thermal 
equilibrium is achieved. 

For the participant matter the compression is much stronger and the excitation 
energy much higher. This provides a short mean free path, but the whole system 
is expanding and cooling fast so that the time available for equilibration is rather 
short and it is questionable if there are enough collisions to ensure local equilibrium 
until freeze out. 



It is very difficult to measure temperature and flow profiles |56|, [57], E8[ and even 



harder or impossible to infer experimentally on the time scale of the evolution of the 
system. Therefore, microscopic transport models which do not assume equilibration 
are needed for a better understanding. These models should go beyond the mean 
field approach, which is a kind of equilibrium assumption in itself, so that in principle 
they are capable to describe many-body correlations like the formation of fragments. 
QMD, AMD and FMD are molecular dynamics models which assert this claim. How 
equilibrium is achieved can then be studied by comparing distributions, for example 
of mass, charge, kinetic energy etc, with equilibrium distributions. 

FMD calculations show that correlations play an important role. Since the time 
evolution of FMD is deterministic, correlations or symmetries can be broken only 
either in the ground state, e.g. in 27 A1 there is no ct-symmetry, but in 28 Si there 
is, or during the time evolution if the Hamilton operator breaks these symmetries 
dynamically. In models like AMD or QMD a randomly fluctuating collision term, 
which models the short-range repulsion as a Langevin force, destroys existing cor- 
relations and symmetries. We do not want to follow this line, because it is unknown 
to which extent correlations of the initial state survive during the collision and it 
even might be an important mechanism for cluster formation in nature. 

In the previous section we investigated beam energies of about 6 AMeV which 
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Figure 10: 19 F- 27 A1 collisions at 32 AMeV, 6 = 0.5 fm: the one-body density in 
coordinate space integrated over z is shown. The contour lines depict the density 
at 0.01, 0.1, 0.5fnT 2 . 



led to dissipative reactions with two out-coming nuclei with about the same mass 
number as in the entrance channel. The energy was not high enough to break neither 
the fused nor the scattered nuclei into pieces. 

Now we chose an energy of 32 AMeV that corresponds to a relative velocity 
between the nuclei of about the Fermi velocity. One expects that at this large 
collective velocity a common mean field cannot be established any more. 
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Figure 11: 19 F- 27 A1 collisions at 32 AMeV, b = 1.5 fm: for explanation see 



fig. 10 



The following figures show a variety of events as contour plots of the one-body 
density in coordinate space. This density is integrated over the ^-direction. Fig- 
ure [11] and [l| present runs at different impact parameters. Within a figure 
different columns show runs with initial states that differ in the orientation of the 



intrinsically deformed ground states. Figure 110 presents four time evolutions at an 
impact parameter of b = 0.5 fm. The time is given in the upper right corners. One 
sees that at this impact parameter and energy a rather long living source is created 
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Figure 12: 19 F- 27 A1 collisions at 32 AMeV, b = 2.5 fm: for explanation see 



which fragments into pieces of all sizes. In fig. |11] for a higher impact parameter 
of b = 1.5 fm the situation is similar but the source is stretched. Figure [T2] finally 
shows runs with an impact parameter of b = 2.5 fm where the angular distribution 
of the fragments is more forward-backward peaked than for the central events where 
it is isotropic. 

Analyzing the time evolution of a cluster one sees that the correlations between 
the wave packets which finally compose the fragment can be followed back for rather 
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long time. 

In the ground state the wave packets arrange in phase space such that the total 
energy is minimized. This FMD ground state is an intrinsic state in which the rela- 
tive positions and orientations of the wave packets reflect many-body correlations. 
If one destroys these correlations by randomly moving all position parameters 7> by 
only 0.2 fm, which leaves the one-body density almost unchanged, a 56 Fe nucleus for 
example achieves already 90 MeV excitation energy. A typical Q-value to break up 
lighter fragments (up to about 16 0) into smaller clusters is 10 MeV. Therefore, for 
not too high excitation energies where the mean free path between hard nucleon- 
nucleon collisions is not yet small compared to the diameter of the fragment, one 
expects the many-body correlations which characterize a bound cluster to survive 
to a large extend. Since the nucleons are indistinguishable one can of course not say 
which nucleon from the projectile or target ends up in a given final fragment. Even 
the wave packets sometimes exchange rapidly r and p, reminding of a Landau-Zener 
crossing. 

It is also interesting to note that in most of the runs several a-particles are 
created (small spheres in the plots), which is in accord with experimental findings 
at a similar energy of E = 35 AMeV in 40 Ca- 40 Ca collisions [[59 . 



At this point a few remarks concerning the time-dependent width might be in 
place, since this degree of freedom is used exclusively in FMD, but not in AMD or 
QMD and led to a discussion about its role [^, |60[1 . 



It was shown not only in our own investigations 03, ESL but also in AMD 



calculations |6T|, |63| , that a time evolution with fixed width parameters does not 
produce evaporation or fragmentation. The reason is that a fixed width implies a 
zero-point energy in each packet of about Skin = 3/ '{AmNCbR) ~ 10 MeV. Inside the 
nucleus this energy is part of the kinetic energy in the Fermi motion. Outside the 
nucleus the zero-point energy of an emitted particle is determined by the amount 
of its localization at the end of the emission process. If the width is kept fixed the 
nucleon carries besides its mean energy always these additional 10 MeV, which is 
much more than the experimentally observed 2 • • • 3 MeV and thus a wave packet 
with fixed width has little probability to escape. 

Compared to that a wave packet with a dynamical width escapes from an excited 
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nucleus by spreading in coordinate space [44], which has two effects: first, the large 



spatial extend reduces the overlap with the other nucleons and the negative potential 
energy tends to zero, second, the positive zero-point energy also becomes small so 
that the sum of both need not change too much during the emission. Classically 
spoken, a particle leaving the nucleus has to climb up the potential wall at the surface 
loosing almost all of its kinetic energy before escaping. FMD gives values of about 
2 MeV for the kinetic energy of evaporated particles in accord with experimental 
findings. 

For fragmentation the dynamical variation of the width parameters is also needed. 
Here a change in the width parameters allows the density to notch so that it may 
break. 

In AMD this lack of evaporation and fragmentation is removed by a subtraction 
of spurious zero-point energy |)2] which also takes care of the localization energy in 
the center of mass motion of the fragments. 



5 Statistical properties of FMD 

Fragmentation reactions show large fluctuations, for example in the mass distribu- 
tion, which are beyond an ensemble averaged mean-field treatment. In this context 
it is important to understand the statistical properties of molecular dynamic models 



especially at low temperatures [p9 | 



There are two aspects. One concerns the thermostatic properties of a molecular 
dynamic model where the attribute thermostatic refers to the properties of the 
static canonical statistical operator, which are contained in the partition function 
Z(T) = Tr(exp{— H/T}). Once the partition function is calculated within a given 
model, its thermostatic properties can be deduced by standard methods like partial 
derivatives of In Z(T) with respect to temperature T or other parameters contained 
in the Hamilton operator H, 

In the case of Fermionic Molecular Dynamics the trace in the partition function 
can be evaluated exactly because the model is based on antisymmetric many-body 
states which form an over-complete set covering the whole Hilbert space. Also the 
states of Antisymmetrized Molecular Dynamics (AMD) [^, |6l| provide a represen- 
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tation for the unit operator. As the calculation of the trace does not depend on 
the representation all thermostatic properties like Fermi-Dirac distribution, specific 
heat, mean energy as a function of temperature etc. ought to be correct and fully 
quantal using FMD or AMD trial states. 

The issue of this section is more the other and even more important aspect, 
namely the dynamical behaviour of a molecular dynamics model. A dissipative 
system which is initially far from equilibrium is expected to equilibrate towards the 
canonical ensemble. The simulation of such a system within the model provides a 
crucial test of its thermodynamic behaviour. 

The time-evolved FMD state is in general not the exact solution of the Schrodinger 
equation, so the correct thermostatic properties do not a priori guarantee correct 
thermodynamic properties. In other words the question is: does the FMD state as 
a function of time explore the Hilbert space according to the canonical weight? 

Since the parameters of the single-particle wave packets follow generalized Hamil- 
ton equations of motion, one is tempted to infer that the dynamical statistical prop- 



erties might be classical |64], ^5|. This conjecture, that classical equations of motion 
always imply classical statistics, is disproved by the following examples, in which 
we compare time-averaged expectation values of wave-packet molecular dynamics 
with the equivalent ones of the canonical ensemble at the same excitation energy. 

Within Fermionic Molecular Dynamics we study the equilibration of four iden- 
tical fermions enclosed in a one dimensional harmonic oscillator. The particles in- 
teract by a weak repulsive two-body potential which is necessary to convert the 
integrable harmonic oscillations into chaotic motion. The important result is that 
the initial many-body state, which is far from equilibrium, approaches the canonical 
ensemble with Fermi-Dirac statistics in an ergodic sense. The time-averaged occu- 
pation numbers of the harmonic oscillator eigenstates are practically identical with 
the Fermi-Dirac distribution of the canonical ensemble, provided the canonical en- 
semble is taken at the time-independent mean excitation energy of the many-body 
state. 

When distinguishable particles, which are described by a product state of wave 
packets, are considered, the molecular dynamic equations for the parameters of the 
wave packets lead to a Boltzmann distribution for the occupation numbers of the 
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single-particle eigenstates. 

A further important result is, that the use of time averages provides us with a tool 
for establishing relations between well-defined quantities of a molecular dynamic 
model such as excitation energy and statistical quantities like temperature. This 



will be done in subsection |5^ where we investigate excited nuclei and the nuclear 
liquid-gas phase transition which is of great experimental 0, [57], Q and theoretical 
interest E3, B7L M M M V% H3 • 



5.1 Thermostatic properties 

In order to describe the thermostatic properties of a given system by means of model 
states exactly, it is necessary that these model states span the whole Hilbert space. 
The question, whether the thermostatic properties of FMD are correct, therefore re- 
duces to the question whether its model states are complete. Starting from coherent 
states the following short explanation proofs that this is the case. 



5.1.1 Completeness relation with coherent states 

Coherent states | z) which are defined as the eigenstates of the harmonic oscillator 
destruction operator a, 

£ I*) = H uo=U j (V a + l) (50) 



form an over-complete set of states. Their completeness relation reads 

(1) _ f d 3 Re (z) d 3 Im (z) \z)(z\ 

J 7T3 (Z\Z) { ) 

d 3 r d 3 p | r, p){f, p\ 



(2tt) 3 (r, p\r, p) 
d 3 r d 3 p \q)(q\ 

WW 

where all three lines are equivalent notations; | z) labeling coherent states by their 
eigenvalue with respect to a, the phase space notation | r, p) labeling the states 
by their expectations values of the coordinate and momentum operators and | q ) 



being the FMD notation. Coherent states are extensively discussed in ref. [74 
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Since we are dealing with fermions the spin degree of freedom has to be considered 
and consequently the resolution of unity changes to 

(1) _ [ d 3 r d 3 p \q)(q\ 
~ ~ J (2*)» ^ (g|g> ' (52) 

where the sum runs over the different magnetic quantum numbers m s which are 
included in the set of parameters denoted by q (see eq. (0)). 

Proceeding one step further the unity operator in the antisymmetric part of the 
two-particle Hilbert space is the antisymmetric product of two single-particle unity 
operators 

1(2) = A (2) (iW® KD) AW (53) 
which may be expressed with antisymmetric two-body states | q±, q2 ) a as 

where | g x , g 2 ) := ^ ( | q x ) ® | q 2 ) - | g 2 ) <8> | gi )) . 

Following this line the resolution of unity in the antisymmetric part of the A-body 
Hilbert space is just the antisymmetric product of single-particle unity operators. 
Be \Q) the unnormalized Slater determinant of single-particle states | q ) and | Q ) 
the normalized Slater determinant 

1 X 

I Q ) = ~X\ sgn ( 7r ) ( I g -(D ) ® " ■ " ® I )) (55) 

7T 

(OIQ>* 

Then the unity operator can be written as the projection of the v4-body unit operator 
onto the antisymmetric subspace of the Hilbert space 

1 {A) = J d/i(Q) \Q){Q\ , (56) 



with a measure 
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that accounts for antisymmetrization by means of the ratio between the norm of 
the Slater determinant (Q\Q) and the norm of the corresponding product state 
Ilfc=i ( Qk | Qk ) (see also J64J). In a sampling where the values of and pu are chosen 
according to Monte Carlo methods this measure determines the probability to find 
the state \Q) belonging to this configuration Q = {fi,p\\ r 2 ,p 2 ; • • • } in Hilbert 
space. If for example two fermions with the same spin are close in r and p then this 
measure is small because the norm ( Q \ Q) = det{{qk\qi)} will be small. 

Coherent states are Gaussian wave-packets with fixed width, but the single- 
particle states of FMD 

f (x — r*) 2 1 
(x\q) oc exp < — \-ip-x>® \x, </>) ® |£) (58) 

contain more degrees of freedom, for instance the width parameter a, the spin angles 
X and and the isospin-3 component £. Since coherent states are already complete 
the additional degrees of freedom a# and ai do not bother. They can be integrated 
keeping track of the normalization as it is shown in (|59|) and the isospin is summed 
over like the spin 

QrQi J (2tt) 3 L R J Qi (q\q) 

where Qr and Qj denote the intervals the width a = a R + iaj is integrated over. 

da R = Q R , / daj = flj . 

n R Jo,! 



The measure then changes to 



dfl(Q) = (Q\Q)[[ — ^ . (60) 



\(qk\qk) (2vr) 3 n R fi/ 



5.1.2 The partition function 



Once the resolution of unity is given in terms of model states the partition function 
can be evaluated. Eq. ([56]) is very useful in calculating traces by means of Monte 
Carlo sampling []64}| . For instance the canonical partition function is given by 

Z(T) = Tr(exp{-£/T}) (61) 
= J dfI(Q) (Q\ exp{-H/T} \Q) . 
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5.1.3 Example 

In the following the above considerations are illustrated with the example of A 



identical fermions in a common one-dimensional harmonic oscillator potential |46 
Starting from the Hamilton operator 

A 



% = E^ H °W' Mn) = ^- + ^mu J 2 x 2 (n) (62) 



n=l 

the mean energy of the v4-fermion system can be derived from the partition function 
Z{T) (fy]) as the derivative with respect to T 

«£»| T = T 2 A ln(Z(T)) (63) 
Jdfl(Q) W(T) Ein O nm {T) [T 2 ^(q m \ exp { - £ho/t} \q n ) 



fdfl(Q)W(T) 

where two the abbreviations W(T) and 0~ X {T) are introduced as 

det((g fc | exp{-/^Ho/T} \qi)\ 
W(T) = — ^ , '-, (64) 

dety{qk\qi)j 

{0-\T)) kl = (q k \ exp{-huo/T} \q t ) . 
The matrix elements can be evaluated in closed form [H6[ . Figure [13| shows the result 



of a Monte-Carlo simulation for a system of four identical fermions in a harmonic 
oscillator with frw = 8 MeV by open circles. The solid line shows the same result, 
but calculating the partition function with eigenstates of the Hamilton operator, 
which obviously gives the same relation. For comparison the classical dependence is 
shown as a dashed-dotted line. 



5.2 Canonical versus ergodic ensemble 

Fermionic Molecular Dynamics is a deterministic microscopic transport theory. Given 
the Hamilton operator and a state | Q(to) ) at a certain time to the state | Q(t) ) is 
known for all times. 

Expectation values are well-defined in FMD so that one can easily calculate 
quantities like the excitation energy of a nucleus or the probability of finding the 
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Figure 13: Excitation energy as a function of temperature, calculated with FMD 
model states (circles) or with eigenfunctions of the harmonic oscillator (solid line). 
The dashed-dotted line shows the classical result. 



system in a given reference state. But it is not obvious how intensive thermody- 
namical quantities, such as the temperature, might be extracted from deterministic 
molecular dynamics with wave packets. In classical mechanics with momentum- 
independent interactions the partition function 

^classical (T) = J JJ - exp j - ^ Classical fal, £l, • • • ) J (65) 

= / n^, x / n<§ -p{-^-->} 

is a product of a term with the kinetic energy and a term containing the interactions. 
Therefore, the momentum distribution can be used to determine the temperature T. 
In the quantum case eq. ( |61~D has to be employed which does not show this factoriza- 
tion. A simple example for this behaviour is the ground state of the free Fermi gas 
where finite momenta are occupied, nevertheless the temperature is zero. Another 
example is the ground state of a nucleus for which the momentum distribution has 
a smeared out Fermi edge due to the finite size and not because of temperature. 

In this section time averaging is compared with a statistical ensemble. If the sys- 
tem is ergodic both are equivalent and statistical properties of molecular dynamics 
can be evaluated by means of time averaging. 
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For this the ergodic ensemble is defined by the statistical operator ^erg as 

R eig := lim — ^— [ * dt \Q(t))(Q(t)\ . (66) 
The ergodic mean of an operator B is given by 
(IT) | { H) ■= Tr (&rg B) = lim _L_ £ dt ( Q(t) | S | Q(t) ) . (67) 

In general the statistical operator Rerg is a functional of the initial state | Q(ti) ), 
the Hamilton operator H and the equations of motion. If the ergodic assumption is 
fulfilled, the statistical operator should only depend on ( H ) , which is actually a 
constant of motion. Thus the average in the ergodic ensemble is always performed 
at the same expectation value of the Hamilton operator. In our notation this is 
denoted by the condition " ( H } " in eq. (|67|) . 



5.2.1 Canonical ensemble of fermions in a harmonic oscillator 

With the statistical operator of the canonical ensemble for A identical fermions (all 
spin up) in a one-dimensional common harmonic oscillator potential ^ho given by 

BT) = J^exp{-%^) (68) 



Z(T) 1 I T 



A oo 



=1 n=0 

the statistical mean of an operator B is calculated as 

« g, »l T -=T (5(T)B) (69) 



|dn^i...d^ (0| B exp (_^E} | 0> 



Z(T) 



E(ni, ■■■ ,n A ) 



Z(T) J ' >'"■»■/ ^ j 

ni<--<nA 

As already mentioned the FMD states are a representation of the unit operator 
and hence can be used to calculate traces. For numerical convenience, however, 
the mathematically identical third line in eq. ( |55[ ) is used, where | ri a , - • - ,n A ) 
denotes the Slater determinant composed of single-particle oscillator eigenstates 
| rii ) , • • • , | n A ) and 

A / 1 \ 

(70) 



£(m, • • • , n A ) = u £ + 2~ ) 
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are the eigenenergies of H no • In eq. ( |69"D the subscript T indicates that the average 
is taken at a constant temperature T. 

In the following a system of four fermions in a common one-dimensional har- 
monic oscillator is investigated. The frequency of the oscillator is chosen to be 
uj = 0.04 fm _1 in order to get a spacing of 8 MeV between the single-particle 
eigenstates. For the canonical ensemble fig. [14] shows the dependence of the excita- 
tion energy on the temperature (l.h.s.) and displays how the lowest eigenstates are 
occupied in the four-fermion system for five different temperatures (r.h.s.). 
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Figure 14: A system of four fermions in a common oscillator described by the 
canonical ensemble. L.h.s.: Excitation energy as a function of temperature (solid 
line). The dashed-dotted line shows the result for a product state (Boltzmann 
statistics). R.h.s.: Occupation numbers p(n) of the oscillator eigenstates for five 
temperatures (eq. (|71|)). The lines are drawn as a guide for the eye. 



The mean occupation probabilities are given by 



(71) 



where c+ denotes the creation operator of a fermion in the oscillator eigenstate 
| n ). 



5.2.2 Ergodic ensemble of fermions in a harmonic oscillator 

In this section the averages of the occupation numbers in the ergodic ensemble 
are evaluated and compared with those of the canonical ensemble discussed in the 
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previous section. As pointed out already, in Fermionic Molecular Dynamics the time 
evolution of Gaussian wave packets in a common oscillator is exact, and thus the 
occupation probabilities of the eigenstates of the Hamilton operator do not change 
in time. In order to equilibrate the system a repulsive short-range interaction Vj is 
introduced. The strength of the interaction is chosen such that the resulting matrix 
elements of Vj are small compared to the level spacing u and the excitation energy 
E*. The contribution of (Vi) to the total energy is of the order of 0.1 . . . 1.0 MeV. 

The initial state is prepared in the following way. Three wave packets with a 
width of a = 1/mio are put close to the origin at x — (—d, 0, d) — with d = 0.5/ ' -Jmw 
- whereas the fourth packet with the same width is pulled away from the centre in 
order to obtain the desired energy. As the mean momenta are all zero, the excitation 
is initially only in potential energy which has to be converted into thermal energy 
by means of the small interaction V\. 

The initial system, which is far from equilibrium, is evolved over about 2000 pe- 
riods of the harmonic oscillator (2tt/uj = 157 fm/c). The equilibration time is rather 
large as we are using a very weak interaction in order not to introduce correlations 
which would destroy the ideal gas picture implied in the canonical ensemble fl68[) of 
non-interacting particles. The time averaging of the occupation numbers (|72|) starts 
at time t\ = 10000 fm/c in order to allow a first equilibration. 



Figure [15] gives an impression of how the occupation numbers evolve in time. 
The part to the left shows the time evolution without interaction which is just a 
unitary transformation in the one-body space. Thus the occupation numbers do not 
change in time although the wave packets are swinging. This has been expected since 
the create eigenstates of the Hamiltonian ^ho- It also serves as an accuracy 
test of the integrating routine. The part to the right displays the evolution with 
interaction at three later times. The occupation probabilities are reshuffled due 
to the interaction and they fluctuate in time. In fig. [16] (l.h.s.) the chaotic time 
dependence of ( Q(t) \ ^n^n \ Q(t) ) for n = 0, 3 and 6 is depicted. 

The result of time averaging is seen in fig. [T7| (symbols) for four different initial 
displacements which correspond to four different excitation energies of the fermion 
system. To each case we assign a canonical ensemble which has the same mean 
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Figure 15: Time evolution of the occupation probabilities for four fermions in a 
common harmonic oscillator potential without (l.h.s.) and with two-body inter- 
action (r.h.s.). The distributions at t — and t = 30000fm/c are connected by a 
solid line. 




t (fm/d xlO 5 n 

Figure 16: L.h.s.: Occupation probabilities versus time — n = 0: circles, n = 3: 
squares, n = 6: triangles. R.h.s.: Variance of the fluctuations Ap 2 (n) calculated 
in the canonical ensemble (solid line) and in the ergodic ensemble (triangles). 



energy. The solid lines in fig. 17 show the corresponding distributions of occupation 
probabilities for these canonical ensembles. Their temperatures T are also quoted 
in the figure. It is surprising to see that there is almost no difference between the 
ergodic and the canonical ensemble: 



(73) 
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Figure 17: Occupation numbers calculated in the ergodic ensemble (symbols, eq. 
([T^)) compared with the canonical ensemble (solid line, eq. (|7T|)). 



provided both have the same excitation energy 



E* = (H UO -E ), {H) = ((H HO -E ))\ T , E = 8u. (74) 



The relation between E* and T is given by eq. (|69|) and displayed in fig. [14 



This result is not trivial because, firstly, the system is very small, consisting 
of only four particles, and secondly, the equations of motion are approximated by 
FMD. The one to one correspondence between the occupation probabilities of the 
ergodic ensemble and the ones of the canonical ensemble, which has the same mean 
energy ( H ) as the pure state, is an impressive demonstration that the system is 
ergodic and that the FMD many-body trajectory covers the phase space according 
to Fermi-Dirac statistics. 



Not only the one-body distributions of the two ensembles coincide, but also the 
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variances of the fluctuations Ap 2 (n), 

Ap 2 H := (( (c+ c n ) 2 »| t - (( c + c n ))| J, , (75) 

as is demonstrated in fig. [16] (r.h.s.). The ergodic mean converges to the result of 
the canonical ensemble which is Ap 2 (n) = p(n)(l — p(n)). 



5.2.3 Canonical and ergodic ensemble for distinguishable particles 

In this section it is shown that time averaging results in quantum Boltzmann statis- 
tics if the fermions are replaced by distinguishable particles. For this end the an- 
tisymmetrized many-body state is replaced by a product state of Gaussian wave 
packets. The resulting equations of motion differ from the FMD case in the skew- 
symmetric matrix A^ U {Q) (given in eq. (p2|) ) which does not couple the generalized 
velocities of different particles any longer. 

For product states the ergodic ensemble is again investigated at different energies 
and compared with the canonical ensembles with the same mean energies. The 
appropriate relation between temperatures and excitation energies in the canonical 
ensemble for distinguishable particles 

E* = ((H no - Eo ))\ T =4| [ooth(^)-l] , E Q = 2u (76) 



is shown by the dashed-dotted line in fig. 14 



Since distinguishable particles are not affected by the Pauli principle, the occupa- 
tion numbers for the many-body ground state look quite different. For instance for 
zero temperature all particles occupy the eigenstate | ) of the harmonic oscillator 



(fig. pi l.h.s 



The initial single-particle states of the interacting system are chosen analogue to 
the fermion case. Again the time evolution of the system exhibits ergodic behaviour 
for all excitation energies. As an example fig. [18] (r.h.s.) is showing the case of 
E*/A = 2.05 MeV (T = 5 MeV) after a time averaging of about 2000 periods. The 
ergodic ensemble (triangles) and the Boltzmann canonical ensemble (solid line) are 
the same within the size of the symbols. The result for Fermi-Dirac statistics with 
the same temperature is included to demonstrate the difference (dotted line). 



Fermionic Molecular Dynamics 



49 




5 10 15 20 5 10 15 20 

n n 



Figure 18: Occupation probabilities for a product state (Boltzmann statistics). 
L.h.s.: Occupation probabilities p(n) of the oscillator eigenstates for five temper- 
atures for the canonical ensemble. R.h.s: Occupation probabilities calculated in 
the ergodic ensemble (symbols) compared with those calculated in the canonical 
ensemble (solid line) for an excitation energy of E* = 2.05 AMeV which corre- 
sponds to a temperature T = 5 MeV in the canonical ensemble. The dotted line 
shows the result for Fermi-Dirac statistics at the same temperature. 

This result shows that equations of motion which are not influenced by the Pauli 
principle lead to the quantum Boltzmann distribution. The only difference to "true 
classical" equations is the presence of the width parameters as dynamical variables. 
Only if they are removed from the equations of motion the statistical behaviour of 
the ergodic ensemble is that of classical statistics. 

5.3 Caloric curve for finite nuclei 

As an outlook of this chapter a method will be presented that allows to investigate 
the caloric curve of finite charged self-bound Fermi systems like nuclei. 

The concept of determining the temperature is to bring a reference system, for 
which thermodynamic relations between temperature and measurable quantities are 
known, into thermal equilibrium with the investigated system. The weakly inter- 
acting ideal gas, where the temperature is given by the mean kinetic energy of the 
particles, may serve as an example. The reference system is called a heat bath if its 
heat capacity is much larger than that of the system and it is called a thermometer 
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if its heat capacity is much less. 




Figure 19: Sketch of the setup: the self-bound excited nuclear system is rep- 
resented by Gaussian wave packets (solid lines) which are enclosed in a broad 
container potential. For the thermometer only one wave packet is shown (dashed 
line), it is contained in a different oscillator. 



As the nuclear system is quantal and strongly interacting its temperature cannot 
be deduced from the momentum distribution or the mean kinetic energy of the 
nucleons. Therefore, the concept of an external thermometer which is coupled to 
the nuclear system is used in the present investigation. The thermometer consists 
of a quantum system of distinguishable particles moving in a common harmonic 
oscillator potential different from the container potential as shown in fig. IT9. 



The time evolution of the whole system is described by the FMD equations of 
motion. For this purpose the many-body trial state is extended and contains now 
both, the nucleonic degrees of freedom and the thermometer degrees of freedom 

\Q}= \Qn)® \Qth) ■ (77) 
The total Hamilton operator including the thermometer is given by 

H = H n {u)+H th + H n _ thl (78) 

where H n {oj) denotes the nuclear Hamiltonian with an additional external field which 
serves as a container. H th is the Hamilton operator of the thermometer system, i.e. 
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the Hamiltonian of a harmonic oscillator. The week repulsive interaction between 
thermometer wave packets and nucleons is given by H n _ th . 

The determination of the caloric curve is done in the following way. The nucleus 
is excited by displacing all wave packets from their ground-state positions randomly. 
Both, centre of mass momentum and angular momentum are kept fixed at zero. To 
allow a first equilibration between the wave packets of the nucleus and those of the 
thermometer the system is evolved over a long time (14000 fm/c). After that a time 
averaging of the energy of the nucleonic system as well as of the thermometer is 
performed over 2000 steps covering a time interval of 2000 fm/c. During this time 
interval the mean of the nucleonic excitation energy 



is evaluated. Eq(N, Z) denotes the FMD ground-state energy of the isotope under 
consideration. The time-averaged energy of the thermometer E t h which is calculated 
during the same time interval determines the temperature T through the relation 
for an ideal gas 



The system is then cooled and after another 2000 fm/c, in which the system equili- 
brates, the averaging is done again. Repeating this procedure one follows the caloric 
curve from high excitations to low excitations. 

The relation between the excitation energy and the temperature is evaluated 
for the three nuclei 16 0, 24 Mg and 40 Ca using the same container potential with 
fiu = 1 MeV. 

The caloric curves shown in the graph on the left hand side of fig. ^U] clearly 
exhibit three different parts. Beginning at small excitation energies the temperature 
rises steeply with increasing energy as expected for the shell model. The nucleons 
remain bound in the excited nucleus which behaves like a liquid drop of fermions. 
At an excitation energy of 3 MeV per nucleon the curve flattens and stays almost 
constant up to about 11 MeV. This plateau at T p& 5 MeV indicates the coexistence 
of liquid and vapour phases, the latter consisting of evaporated nucleons which 
are in equilibrium with the residual liquid drop due to the containment. Around 



E* 




(79) 



T 




(80) 
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Figure 20: L.h.s.: caloric curve of 16 0, 24 Mg and 40 Ca at ftw — 1 MeV, r.h.s.: 
time averaged radial density distribution of 24 Mg at various excitation energies 
in the coexistence region. 



E* / A 11 MeV all nucleons are unbound and the system has reached the vapour 
phase. This is indicated by the steep rise of the caloric curve beyond this point. 
One has to keep in mind that the plateau is not the result of a Maxwell construction 
as in nuclear matter calculations. In the excitation energy range between 3 and 
11 MeV per particle an increasing number of nucleons is found in the vapour phase 
outside the liquid phase. This can be seen in the density plot on the right hand 
side of fig. ^ where the radial dependence of the time-averaged density for 24 Mg 
at three excitation energies in the coexistence region is shown. For small excitations 
(E*/A = 4.0 MeV) the nucleus, which due to recoil is bouncing around, is surrounded 
by very low density vapour (solid line). The dashed line (E*/A = 7.4 MeV) and 
the dashed-dotted line (E*/A = 11.2 MeV) show that with increasing energy the 
vapour contribution is growing and the amount of liquid decreasing. However, in 
the high energy part of the plateau the averaged one-body density shown here does 
not represent the physical situation adequately. The time-dependent many-body 
state shows the formation and disintegration of several small drops, which due to 



time-averaging cannot be seen in fig. ^0[ Above E*/A ~ 13 MeV only vapour is 
observed. 

In how far the experimental result, fig. [21], where the temperature is determined 
from isotope ratios |[56||, can be compared to our time-averaged temperature and 
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Figure 21: Caloric curve extracted from spectator fragmentation by the ALADIN 



collaboration (picture taken from ref. [56]) 



excitation energy is still under investigation. 

The results on the caloric curve of nuclei present another example for the wide 
applicability of Fermionic Molecular Dynamics. A correct equilibrium behaviour is 
of course a necessary condition for further investigations of non-equilibrium situa- 
tions as in nucleus-nucleus collisions. These investigations are followed up presently. 
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A Interaction matrix elements 



A.l Approximation of the matrix elements 



A large fraction of the numerical effort of FMD is caused by the evaluation of the 
interaction matrix elements 



V(Q)= (V) =Tr( Z P^) 



1 A 

k,l,m,n 



(ikQil £ \q m q n ) { O mk O n i - O m iO. 



nk 



(81) 



This effort grows with A 4 because Gaussian wave packets are not orthogonal. In 
order to reduce the computation time two types of approximations were tested [45. 
46| . Both approximation schemes use that antisymmetrization effectively reduces 
the strength of the interaction. 



The first ansatz [45 



V(Q) 



(?fc?l| Z \lkqi)a 

y — / — \ — \ — ex P 



k<l 



cm 



(qkQMkQl) 



"4 C « 



(82) 



E 



(qkq m \qmqk) (qiq m \q m qi) 
(qkq m \qkq m ) (qiq m \qiq m } 



takes all matrix elements (q k qi\ V \qkQi) an d {q k qi\ V \qiqk) into account which are 
scaled with an overall factor counting the overlapping wave packets. 



The second ansatz j|6| reduces this effort even further since it needs only the 
direct matrix elements 



V(Q) 



in 



f£ {qkqMkqi) 



(1 - Iki) + 7fc/exp 



(qkqMiqk} 
(qkqMkqi) 



(83) 



erf { - 

2 



A 

E 

mj^k,l 



(qkqm\q m qk) (qiqm\q m qi) 
(qkq m \qkq m ) (qiq m \qiq m } 



Both schemes (^) and (|83| ) provide a reasonable approximation for the expectation 
value of the two-body interaction. In order to reproduce the nuclear binding energies 
it can be necessary to readjust the strength of the original interaction. 
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A. 2 Coulomb interaction 

This section briefly shows how the second approximation fl8"3|) works for the Coulomb 



interaction (f40|) 



i ^ —* ir / • . \ I t _ . 1.44MeVfm p n PrV^ ^uv-. -\ /o,x 

(ar^Xj-l £ = -p ^t-P ®P 8\x i -x k )5\x 1 -x l ) (84) 



| *^ 2 J I 



where P P denotes the projection operator on the protons. The diagonal matrix 



elements are 



' MW ' Vfm -feWiS (85) 



with 



(qkqMkqi) r kl 



r k i = \r k - n | (86) 
2 ct/cij am 



a fejR |a z | 2 + a m |a fc | 2 



2 P 

erf(x) := — p= / du exp{— u 2 } 
V 71 " Jo 



Pli = (&l£ p l6) = i(i + 6X&l6> (87) 

The isospin variable takes values ^ = 1 for protons and = —1 for neutrons, 
respectively. 

Figure [E] shows the result of the approximation for two systems, two identical 
protons and a 12 C nucleus. The approximation (thin line) compares nicely to the 
exact result (thick line). 
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Figure 22: Approximation of the Coulomb energy. The l.h.s shows the Coulomb 
energy of two protons with the same spin component as a function of the relative 
distance, the r.h.s. shows the Coulomb energy of a 12 C nucleus as a function 
of the relative distance between the three alpha clusters. The thick solid line 
represents the exact solution, the dashed line the direct term only and the thin 
line the result of the approximation. 
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